A Nonconservative Earthquake Model of Self-Organized Criticality on a Random 

Graph 
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We numerically investigate the Olami-Feder-Christensen model on a quenched random graph. 
Contrary to the case of annealed random neighbors, we find that the quenched model exhibits 
self-organized criticality deep within the nonconservative regime. The probability distribution for 
avalanche size obeys finite size scaling, with universal critical exponents. In addition, a power law 
relation between the size and the duration of an avalanche exists. We propose that this may represent 
the correct mean-field limit of the model rather than the annealed random neighbor version. 

PACS numbers: 05.65. +b, 45.70.Ht, 89.75.-k 

The idea of self-organized criticality (SOC) was intro- 
duced as a possible explanation for the widespread oc- 
currence in nature of long range correlations in space 
and time Jl]]. The term refers to the intrinsic tendency 
of a large class of spatially extended dynamical systems 
to spontaneously organize into a dynamical critical state. 
In general, SOC systems are driven externally at a very 
slow rate and relax with bursts of activity, avalanches, 
on a very fast time scale. One signature of SOC is a 
scale free, e.g. power law, distribution of avalanche sizes. 
This is normally related to some long range spatial and 
temporal correlations within the system. Typical natu- 
ral realizations of this phenomena include, among others, 
earthquakes, forest fires, and biological evolution (for re- 
views, see (2^). 

A problem that has attracted a lot of attention, but 
is still poorly understood, is that of identifying funda- 
mental mechanisms leading to SOC behavior. In partic- 
ular, much effort has been directed at understanding how 
conservation of the transported quantity (e.g. sand) in 
the avalanche dynamics affects criticality ^|||. For in- 
stance, it is well known that the Abelian sandpile model 
ijlj, is subcritical if dissipation is introduced ||. On 
the other hand, nonconservative sandpile models which 
display criticality have since been introduced Al- 
though both the analytical and numerical evidence in 
favor of criticality are quite convincing, the role played 
in these models by the non-conservative dynamics is not 
clear. In fact, dissipation is a dynamical variable and 
does not always occur. 

A model, which in the context of SOC in nonconserva- 
tive systems has played an important role, is the Olami- 
Feder-Christensen (OFC) model of earthquakes M. In 
the OFC model a finite fraction, controlled by a fixed 
parameter a, of the transported quantity is dissipated in 
each relaxation event. The presence of criticality in the 
non-conservative OFC model has been controversial since 
the very introduction of the model |l^] and it is still de- 
bated 0,0. Recent numerical investigations, though, 
have shown that the OFC model on a square lattice dis- 
plays scaling behavior, up to lattice sizes presently ac- 



cessible by computer simulations [ |13|Jl4| ] . The avalanche 
size distribution is described by a power law, character- 
ized by a universal exponent t ~ 1.8, independent of the 
dissipation parameter. This distribution does not display 
finite size scaling, however. 

To overcome the limitation of relying almost exclu- 
sively on computer simulation results, it has sometimes 
been useful to consider an annealed random neighbor 
(RN) version of the model JToj-p^l, where each site in- 
teracts with randomly chosen sites instead of its nearest 
neighbors on the lattice. This considerably simplifies the 
problem. In the past, RN models have usually been con- 
sidered as mean-field descriptions of their fixed lattice 
counterparts, since spatial correlations are absent. Anal- 
ogous to other RN models, the RN OFC model can be 
solved analytically |l6|,[l7]]. It displays criticality only in 
the conservative case, where it becomes equivalent to a 
critical branching process. As soon as some dissipation 
is introduced the avalanches become localized, although 
the mean avalanche size diverges exponentially as dissi- 
pation tends to zero. The absence of criticality, together 
with the exponential divergence has cast some doubt on 
whether the OFC model on a fixed lattice is critical. 

However, it is important to point out that the RN 
model may not describe the behavior of the OFC model 
on a fixed lattice in any dimension, and thus may not 
correspond to the mean field limit of the model. Usually, 
mean field behavior describes the high dimensional be- 
havior of the system (e.g. the behavior above an upper 
critical dimension); this is not exactly the same limit as 
a model without any spatial correlations in it. 

In fact, criticality in the OFC model on a lattice has 
been ascribed to a mechanism of partial synchroniza- 
tion The system has a tendency to order into a 
periodic state |l9|-[2l|] which is frustrated by the presence 
of inhomogeneities such as the boundaries. In addition, 
inhomogeneities induce partial synchronization of the el- 
ements of the system building up long range spatial cor- 
relations and thereby creating a critical state. The mech- 
anism of synchronization requires an underlying spatial 
structure and therefore cannot operate in an annealed 



1 



RN model, where each site is assigned new random neigh- 
bors at each update. 

The main purpose of the present work is the investi- 
gation of the OFC model on a quenched random graph. 
This formulation, which can be handled numerically, is 
worth analyzing to see if it presents critical or non-critical 
behavior. Since the largest distance between two sites in 
the random graph grows only as a logarithm of the num- 
ber of sites, it can be considered to be a high dimensional 
limit of a lattice model, and thus may describe the cor- 
rect mean held limit. Contrary to the RN case, in a ran- 
dom graph the choice of neighbors is not annealed but 
quenched, so that spatial correlations can develop. In- 
deed, we show that the OFC model on a random graph 
displays criticality even in the nonconservative regime. 

A random graph is defined as a set of N sites connected 
at random by bonds. Two connected sites are denoted 
as "nearest-neighbor" . Formally, the random graph can 
be constructed by considering all N(N — l)/2 possible 
bonds between sites and occupying a certain number of 
them with equal probability. A constraint of fixed con- 
nectivity can also be imposed by requiring that each site 
has the same number of neighbors, q. We have mainly 
concentrated on this latter situation but the first case will 
also be discussed. The model is then defined as follows. 
To each site of the graph is associated a real variable Fi, 
which initially takes some random value in the interval 
(0, Fth)- All the forces are increased uniformly and simul- 
taneously at the same speed, until one of them reaches 
the threshold value F t h and becomes unstable (Fi > F t h)- 
The uniform driving is then stopped and an "earthquake" 
(or avalanche) starts: 



Fi > F th 



Fi^O 



aF 



(1) 



where "nn" denotes the set of nearest-neighbor sites of 
i. The parameter a controls the level of conservation of 
the dynamics and, in the case of a graph with fixed con- 
nectivity q, it takes values between and 1/q (a = 1/q 
corresponding to the conservative case). The toppling 
rule (Q) can possibly create new unstable sites, produc- 
ing a chain reaction. All sites that are above threshold at 
a given time step in the avalanche relax simultaneously 
according to (|l|) and the earthquake is over when there 
are no more unstable sites in the system (Fi < F t h, Vi). 
The uniform growth then starts again. The number of 
topplings during an earthquake defines its size, s, and we 
will be interested in the probability distribution P/v(s). 
Another quantity of interest is the duration t of an earth- 
quake which will be identified with the number of time 
steps needed for the earthquake to finish. 

We consider first a random graph where all sites have 
exactly the same number of nearest neighbors q. In this 
case, we have verified (both for q = 4 and q = 6) that 
the system organizes into a subcritical state. This is 



analogous to what happens in the OFC model on a lat- 
tice with periodic boundary conditions, where no critical 
behavior is observed [jl9-21|. In order to observe scal- 
ing in the avalanche distribution, one has to introduce 
some inhomogeneities. In the lattice model this is gen- 
erally achieved by considering open boundary conditions 
which imply that boundary sites have fewer neighbors 
and therefore cycle at a different frequency from bulk 
sites. This is an inhomogeneity with a diverging length 
scale in the thermodynamic limit. For the OFC model on 
a random graph, we have found that it suffices to consider 
just two sites in the system with coordination q — 1 [p2l. 
When either of these sites topple according to rule (|]), 
an extra amount aFi is simply lost by the system. 

After a sufficiently long transient time, the system set- 
tles into a statistically stationary state. We have ver- 
ified that the statistical properties of the system (e.g. 
the avalanche distribution) are independent of the actual 
realization of the random graph, as long as the coordi- 
nation number q is the same. As a point of compari- 
son, in figure 1 we report the probability distribution of 
avalanche sizes for (a) the annealed RN model and (b) 
for the OFC model on a random graph for various sys- 
tem sizes N. The dynamical rule for the annealed RN 
model are formally similar to (1), where, instead of the 
nearest-neighbor sites on the graph, q new random sites 
are chosen at each relaxing event. In both cases of fig. 1, 
the number of neighbors is q = 4 and the parameter 
a = 0.10. It is clear that no scaling is present in the 
RN model as the cut-off in the avalanche size distribu- 
tion does not grow with system size. On the contrary 
for the model on a random graph, the distribution scales 
with system size, which is indicative of a critical state. In 
fact, the largest avalanche roughly coincides with system 
size. It is important to underline that we are considering 
a situation far away from the conservative case (60% of 
the force in the toppling site is dissipated) and therefore 
one could expect that if a finite length scale related to 
conservation existed in the system it should appear for 
system sizes we have considered. 

In order to characterize the critical behavior of the 
model, a finite size scaling (FSS) ansatz is used, i.e. 

(2) 



P N (s) ~ N-Pf{s/N D ) 



where / is a suitable scaling function and f3 and D are 
critical exponents describing the scaling of the distribu- 
tion function. In figure 2, a FSS collapse of P/v(s) for 
different values of a and for different q is shown. The 
distribution Pn{s) satisfies the FSS hypothesis reason- 
ably well, with universal critical coefficients. The criti- 
cal exponent derived from the fit of fig. 2 are (3 ~ 1.65 
and D = 1, independent of the dissipation parameter 
a and the coordination number of the graph q. The 
FSS hypothesis implies that, for asymptotically large 
N, Pn(s) ~ s~ T and the value of the exponent is 
t = /3/D ~ 1.65. Due to the numerical uncertainty on 
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the estimate it is difficult to assert with certainty that r 
is a novel exponent, different from the one for the con- 
servative RN model (r = 1.5) or the lattice model in two 
dimensions (r ~ 1.8). 

The OFC model on a lattice does not show ordinary 
FSS Although the avalanche size distribution con- 
verges to a well-defined, universal power law, the cut-off 
in the distribution due to finite system sizes does not be- 
have according to FSS |l3). In particular, the apparent 
numerical value for the exponent D determined through 
FSS would violate some exact bounds In fact, in the 
non-conservative model each site can only discharge a fi- 
nite number of times during an avalanche, which imposes 
D < 1 in eq. (0) . In order to recover standard FSS in the 
two dimensional lattice model one has to consider earth- 
quakes localized within subsystems of linear dimension 
small compared to the overall system size fl4|| . 

We have performed numerical simulations of the OFC 
model on a square two dimensional periodic lattice where 
just two sites have 3 neighbors. We have verified that the 
avalanche size distribution scales with system size but, as 
for the OFC model with open boundary conditions, FSS 
appears to be violated in the cut-off region. The power 
law exponent of the distribution is consistent with the 
exponent of the OFC model with open boundary con- 
ditions, i.e. t ~ 1.8. For the range of system sizes we 
could simulate, the critical behavior of the model on a 
lattice and on a random graph (with the same number 
of defects) appear to be different. In particular SOC on 
the quenched random graph appears to be described by 
ordinary FSS, consistent with a mean field limit. 

We now discuss the time properties of the avalanches. 
In fig. 3 we report the average size of an avalanche stop- 
ping at time t, < s >t, as a function of the rescaled time 
t = t + 2 (as we are mainly interested at large values of t, 
the constant should be irrelevant). The curves for differ- 
ent system sizes overlap (deviations can be attributed to 
finite-size effects) and we observe that < s >t— t 1 , where 
7 ~ 2.1, providing further evidence of criticality in the 
nonconservative system. 

An interesting question (of difficult solution though) is 
whether the model becomes subcritical below a certain 
non-zero value a c (for a = the system is clearly not 
critical as sites do not interact). In our simulations we 
have found that for a > 0.10 and q = 4 the model dis- 
plays scaling behavior with universal critical exponents 
(sec fig. 1 and 2). For lower values of a the analysis is 
more complicated as the extremely long transient times 
to stationarity prevent the investigation of large lattices. 
Nonetheless for very low values of a (a ~ 0.03) the cut- 
off in the avalanche distribution does not appear to vary 
systematically with system size (even for relatively small 
systems), suggesting that a non-zero a c may exist. 

We have also considered the OFC model on a random 
graph with variable local connectivity qi. In this case, the 
toppling rule (fw must be modified to take into account 



that different sites have a different coordination number 
qi. Each site consequently has a different oti, which we 
determined by requiring that the total fraction a of the 
force transferred from the unstable site to the nearest- 
neighbors sites is constant in the system, i.e. = a/qt. 
In particular, we have studied a graph with average con- 
nectivity < qi >= 4. We have found that there is no crit- 
icality in the system as the cut-off in the probability dis- 
tribution does not scale with system size. In agreement 
with previous investigations p3| , p4[ , this result indicates 
that if the disorder is too strong (as for a completely 
random graph) the critical state is destroyed. On the 
other hand inhomogeneities are necessary to break the 
periodic state the system would otherwise reach. It is a 
difficult question to establish what is the maximum level 
of disorder that the system can sustain without loosing 
its critical properties. 

In conclusion, in this paper we have investigated the 
OFC model on a quenched random graph. We have 
shown that the model is critical even in the nonconser- 
vative regime. This is in contrast to what happens in 
the annealed RN OFC model which displays criticality 
only in the conservative case. Contrary to the annealed 
case, a quenched random graph has an underlying spatial 
structure so that partial synchronization of the elements 
of the system can still occur. As a random graph can be 
regarded as a high dimensional limit of a regular lattice, 
we propose that this formulation represents the correct 
mean-field limit of the model rather than the annealed 
random neighbor version. 
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FIG. 2. Finite-size scaling plots for Pjv(s) for (a) q = 4, 
a = 0.15, (b) q = 4, a = 0.20 and (c) q = 6, a = 0.10. System 
sizes are N = 4- 10 3 , 16- 10 3 , 64- 10 3 and 256- 10 3 . The critical 
exponents are (3 = 1.65 and D = 1. For visual clarity, curves 
(a) and (c) have been shifted along the x axis, x — ► x — 1 and 
x — » x + 1, respectively. 
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FIG. 1. Probability distribution (a) for the RN OFC 
model and (b) for the OFC model on a random graph. In 
both cases, q = 4 and a = 0.10. System sizes are (a) N = 10 3 , 
4 • 10 3 and (b) N = 10 3 , 16 • 10 3 , 256 • 10 3 . 
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FIG. 3. Average size of an avalanche lasting t time steps 
as a function oft for q — 4 and (a) a = 0.15 and (b) a = 0.20. 
Different curves correspond, from bottom to top, to system 
sizes N = 1 • 10 3 , 4 • 10 3 , 16 • 10 3 , 64 • 10 3 and 256 • 10 3 . 
The slope of the straight line is 7 = 2.1. Curve (b) has been 
shifted along the x axis, x = x + 1, for visual clarity. 
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